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Abstract 

Many questions in physical cosmology regarding the thermal history of the intergalactic 
medium, chemical enrichment, reionization, etc. are thought to be intimately related to the 
nature and evolution of pregalactic structure. In particular the efficiency of primordial star for- 
mation and the primordial IMF are of special interest. We present results from high resolution 
three-dimensional adaptive mesh refinement simulations that follow the collapse of primordial 
molecular clouds and their subsequent fragmentation within a cosmologically representative vol- 
ume. Comoving scales from 128 kpc down to 1 pc are followed accurately. Dark matter dy- 
namics, hydrodynamics and all relevant chemical and radiative processes (cooling) are followed 
self-consistently for a cluster normalized CDM structure formation model. Primordial molecular 
clouds with ~ 10^ solar masses are assembled by mergers of multiple objects that have formed 
hydrogen molecules in the gas phase with a fractional abundance of < 10~^. As the subclumps 
merge cooling lowers the temperature to ~ 200 K in a "cold pocket" at the center of the halo. 
Within this cold pocket, a quasi-hydrostatically contracting core with mass ~ 200 Mq and num- 
ber densities >, 10^ cm"'^ is found. We find that less than 1% of the primordial gas in such small 
scale structures cools and collapses to sufficiently high densities to be available for primordial star 
formation. Furthermore, it is worthwhile to note that this study achieved the highest dynamic 
range covered by structured adaptive mesh techniques in cosmological hydrodynamics to date. 



1. Introduction 

Saslaw and Zipoy (1967) realized the impor- 
tance of gas phase H2 molecule formation in pri- 
mordial gas for the formation of proto-galactic ob- 
jects. Employing this mechanism in Jeans unsta- 
ble clouds, Peebles and Dicke (1968) formulated 
their model for the formation of primordial glob- 
ular clusters. Further pioneering studies in this 
subject were carried out by Takeda etal. (1969), 
Matsuda etal. (1969), and Hirasawa etal. (1969) 
who followed in detail the gas kinetics in collapsing 
objects and studied the possible formation of very 
massive objects (VMO's). In the 1980's the pos- 
sible cosmological consequences of population III 



star formation were assessed (Rees and Kashlinsky 
1983; Carr etal. 1984; Couchman and Rees 1986). 
In particular Couchman and Rees (1986) discussed 
first structure formation within the standard cold 
dark matter model. Their main conclusions were 
that the first objects might reheat and reionize the 
universe, raise the Jeans mass and thereby influ- 
ence subsequent structure formation. 

Early studies focused on the chemical evolu- 
tion and cooling of primordial clouds by solving 
a chemical reaction network within highly ideal- 
ized collapse models (cf. Hirasawa 1969; Hutchins 
1976; Palla etal. 1983; MacLow and ShuU 1986; 
Puy etal. 1996; Tegmark etal. 1997). Some hy- 
drodynamic aspects of the problem were stud- 



ied in spherical symmetry by Bodenheimer (1986) 
and Haiman, Thoul and Loeb (1996). Recently 
multi-dimensional studies of first structure forma- 
tion have become computationally feasible (Abel 
1995; Anninos & Norman 1996; Zhang etal. 1997; 
Gnedin & Ostriker 1997; Abel etal. 1998a, 1998b; 
Bromm etal. 1999). These investigations have 
provided new insights into the inherently mul- 
tidimensional, nonlinear, nonequilibrium physics 
which determine the collapse and fragmentation 
of gravitationally and thermally unstable primor- 
dial gas clouds. 

In Abel, Anninos, Norman & Zhang 1998a 
(hereafter AANZ) we presented the first self- 
consistent 3D cosmological hydrodynamical sim- 
ulations of first structure formation in a stan- 
dard cold dark matter-dominated (SCDM) uni- 
verse. These simulations included a careful 
treatment of the formation and destruction of 
H2 — the dominant coolant in low mass halos 
(Mtot = 10^ - 10* Mq) which collapse at high 
redshifts (z ^ 30 — 50). Among the principal find- 
ings of that study were: (1) appreciable cooling 
only occurs in the cores of the high density spher- 
ical knots located at the intersection of filaments; 
(2) good agreement was foimd with semi-analytic 
predictions (Abel 1995; Tcgmark etal. 1997) of the 
minimum halo mass able to cool and collapse to 
higher densities; (3) only a small fraction (< 10%) 
of the bound baryons are able to cool promptly, 
implying that primordial Pop III star clusters may 
have very low mass. Due to the limited spatial 
resolution of the those simulations 1 kpc co- 
moving), we were unable to study the collapse to 
stellar densities and address the nature of the first 
objects formed. 

In this paper we present new, higher- resolution 
results using the powerful numerical technique of 
adaptive mesh refinement (AMR, Bryan & Nor- 
man 1997; Norman & Bryan 1999) which has shed 
some light on how the cooling gas fragments. With 
an effective dynamic range of 262, 144 the numeri- 
cal simulations presented here are the highest reso- 
lution simulations in cosmological hydrodynamics 
to date. Although we are not yet able to form 
individual protostars, we are able to resolve the 
collapsing protostellar cloud cores which must in- 
evitably form them. We find the cores have typi- 
cal masses ^ 2OOM0, sizes ~ 0.3 pc, and number 
densities n > 10^ cm~'^-similar to dense molec- 



ular cloud cores in the Milky Way with one vi- 
tal difference: the molecular hydrodgen fraction is 
~ 5 X 10~*, meaning the cores evolves very differ- 
ently from Galactic cores. 

The plan of this paper is as follows. The simu- 
lations are briefly described in Sec. 2. Results are 
presented in Sec. 3. The properties and fate of 
the primordial protostellar cloud are discussed in 
Sec. 4. Conclusions follow in Sec. 5. Results of a 
broader survey of simulations will be reported in 
Abel, Bryan & Norman (1999). 

2. Simulations 

The three dimensional adaptive mesh refine- 
ment calculations presented here use for the hydro- 
dynamic portion an algorithm very similar to the 
one described by Berger and Collela (1989). The 
code iitilizes an adaptive hierarchy of grid patches 
at various levels of resolution. Each rectangular 
grid patch covers some region of space in its parent 
grid needing higher resolution, and may itself be- 
come the parent grid to an even higher resolution 
child grid. Our general implementation of AMR 
places no restriction on the number of grids at a 
given level of refinement, or the number of levels 
of refinement. However, we do restrict the refine- 
ment factor - the ratio of parent to child mesh 
spacing - to be an integer (chosen to be 2 in this 
work). The dark matter is followed with meth- 
ods similar to the ones presented by Couchman 
(1991). Furthermore, the algorithm of Anninos 
etal. (1997) is used to solve the time-dependent 
chemistry and cooling equations for primordial gas 
given in Abel etal. (1997). More detailed descrip- 
tions of the code are given in Bryan & Norman 
(1997, 1999), and Norman & Bryan (1999). 

The simulations are initialized at redshift 100 
with density perturbations of a SCDM model with 
flB = 0.06, h = 0.5, and erg = 0.7. The abun- 
dances of the 9 chemical species (H, H+ , H~ , He, 
He+ , He++ , H2 , Hj" , e~) and the temperature 
are initialized as discussed in Anninos and Nor- 
man (1996). After a collapsing high-cr peaks has 
been identified in a low resolution run, the simula- 
tion is reinitialized with multiple refinement levels 
covering the Langrangian volume of the collaps- 
ing structure. The mass resolution in the initial 
conditions within this region are 0.53(8.96) Mq in 
the gas (dark matter) . The refinement criteria en- 
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sure that: (1) the local Jeans length is resolved 
by at least 4 grid zones, and (2) that no cell con- 
tains more than 4 times the initial mass element 
(0.53 Mq). Wc limit the refinement to 12 levels 
within a 64"^ top grid which translates to a maxi- 
mum dynamic range of 64 x 2^^ = 262, 144. 

3. Results 

We find that primordial molecular clouds are 
only formed at the intersection of filaments, in 
agreement with the results of AANZ. The evo- 
lution of these primordial molecular clouds is 
marked by frequent mergers yielding highly com- 
plex velocity and density fields within the "virial" 
radius. In the following three sections, we first 
describe the evolution of these objects and then 
their morphology and structure. 

3.1. Formation of the First Objects 

To illustrate the physical mechanisms at work 
during the formation of the first cosmological ob- 
ject in our simulation, we show the evolution of 
various quantities in Figure 1. The top panel of 
this plot shows the virial mass of the largest object 
in the simulation volume. We divide the evolution 
up into four intervals. In the first, before a redshift 
of about 35, the Jeans mass in the baryonic com- 
ponent is larger than the mass of any non-linear 
perturbation. Therefore, the only collapsed ob- 
jects are dark-matter dominanted, and the bary- 
onic field is quite smooth. (We remind the reader 
that a change in the adopted cosmological model 
would modify the timing, but not the nature, of 
the collapse.) 

In the second epoch, 23 < z < 35, as the non- 
linear mass increases, the first baryonic objects 
collapse. However, these cannot efficiently cool 
and the primordial entropy of the gas prevents 
dense cores from forming. This is shown in the 
second frame of figure 1 by a large gap between the 
central baryonic and dark matter densities (note 
that while the dark matter density is limited by 
resolution, the baryonic is not, so the true difi^er- 
ence is even larger). As mergers continue and the 
mass of the largest clump increases, its temper- 
ature also grows, as shown in the third panel of 
this figure. The H2 fraction also increases (bot- 
tom panel). 

By 2; 23, enough H2 has formed (a few 



xlO~^), and the temperature has grown suffi- 
ciently high that cooling begins to be important. 
During this third phase, the central temperature 
decreases and the gas density increases. How- 
ever, the collapse is somewhat protracted because 
around this point in the evolution, the central den- 
sity reaches n ^ 10'* cm^^, and the excited states 
of H2 are in LTE. This results in a cooling time 
which is nearly independent of density rather than 
in the low density limit where tcool ^ (e.g. 
Lepp & ShuU 1983). 

Finally, at 2; ~ 19, a very small dense core forms 
and reaches the highest resolution that we allowed 
the code to produce. It is important to note that 
at this point, the maximum gas density in the sim- 
ulations exceeds 10* cm""*, and at these densities, 
3-body formation of molecular hydrogen will be- 
come dominant (see Palla etal. 1983). Also, the 
assumption of optical thin cooling begins to break 
down and radiative transfer effects become impor- 
tant. Therefore, only simulation results at and 
above this redshift will be discussed. It is worth- 
while to note that the simulations presented here 
are physics rather than resolution limited. 

3.2. Morphology 

The increase in dynamic range by ~ 1000 in the 
simulations presented here as compared to AANZ 
allow us to investigate the details of the fragmen- 
tation process in detail. Visualizations of the gas 
density and temperature on two different scales at 
z = 19.1 are shown in the color plate 3. In the 
upper left panel the velocity field is shown super- 
imposed on the density. The 5 x 10^ Mq struc- 
ture forms at the intersection of two filaments with 
overdensities of ~ 10. Most of the mass accretion 
occurs along these filaments. The complexity of 
the velocity field is evident; the accretion shock is 
highly aspherical and of varying strength. Within 
the virial radius (r = 106 pc), there are a number 
of other cooling regions. The right-hand panels 
zooms in on the collapsing fragment, (note that 
the smallest resolution element (0.02 pc) in the 
simulations is still 1600 times smaller than the slice 
shown in the right panels). The small fragment in 
the center of this image has a typical overdensity 
of > 10*^ and a mass of ~ 200 M„. 
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Fig. 1. — The top panel shows the evolution of the virial mass of the most massive clump as a function 
of redshift. The remaining panels show the density (both dark and baryonic), the temperature, and the 
molecular hydrogen mass fraction at the central point of that clump. The central point is defined as the 
point with the highest baryon density. Clearly the finite gas pressure prevents baryons from clumping as 
much as the dark matter at redshifts >, 23. The vertical line at z = 19.1 indicates where our numerical 
model breaks down. 
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3.3. Profiles 

Despite the complex structure of tlie primor- 
dial molecular clouds much of their structure can 
be understood from spherical profiles of the phys- 
ical quantities, particularly for the dense central 
core which is nearly spherical. Figure 3.3 shows 
mass- weighted, spherical averages of various quan- 
tities around the densest cell found in the simula- 
tion at redshift 19.1. Panel a) plots the baryon 
number density, enclosed baryon mass, and local 
Bonnor-Ebert mass-' w 27 MqT}^^ / ^Jn versus ra- 
dius. Panel b) plots the abundances of H2 and free 
electrons. Panel c) compares three timescales de- 
fined locally: the H2 cooling time tn^ , the freefall 
time tff = [37r/(32Gp)]^/^, and the sound cross- 
ing time tcross = r/cs = 7.6 x lQ^rpc/\/TK yrs. In 
panel d) we identify two distinct regions — labeled 
I and II — as defined by the temperture profile. 
Region I ranges from outside the virial radius to 
rT„.„ ^ 5 pc, the radius at which the infalling 
material has cooled down to 

Tmin ^ 200 K— near 
the minimum temperature allowed by H2 cooling. 
Within region I, the temperature profile reflects, 
in order of decreasing radius, cosmic infall, shock 
virialization, adiabatic heating in a settling zone, 
and an H2 cooling flow. In region II, the temper- 
ature slowly rises from Tmin to ~ due to 
adiabatic heating. 

For most of region I the II2 cooling time is 
comparable to the free-fall time, as is illustrated 
in panel c) of Figure 3.3. The H2 number fraction 
rises from 7 x 10~^ to 2 x 10^'* as the free elec- 
tron fraction drops from 2 x 10~^ to 2 x 10~^. At 
rT„i„, the sound crossing time becomes substan- 
tially shorter than the cooling time. This suggests 
that region II is contracting quasi-hydrostatically 
on the cooling time scale, which approaches its 
constant high-density value at small radii. This 
constant cooling time of ~ 10^ years sets the time 
scale of the evolution of the fragment until it can 
turn fully molecular via three body associations. 
Inside r 0.3 pc, the enclosed baryonic mass of 
^ 2OOM0 exceeds the local Bonnor-Ebert mass, 
implying this material is gravitationally unstable. 
However, due to the inefficient cooling, its collapse 
is subsonic (panel e). The radius where M > Mbe 



Bonnor-Ebert mass is the analog of the Jeans Mass but 
assuming an isotlicrmal {p oc r~'^) instead of a uniform 
density distribution. 



defines our protostellar cloud core. 

4. Discussion 

Many interesting features of the collapsing and 
fragmenting "primordial molecular cloud" are 
identified. Most notable is the formation of an 
initially quasi-hydrostatically contracting core of 
~ 200 Mq which becomes gravitationally unsta- 
ble. We argue that this is a characteristic mass 
scale for core formation mediated by H2 cool- 
ing. Substituting into the formula for the Bonnor- 
Ebert mass Tmin and tllte we get 240Mq. 

What will be the fate of the collapsing core? 
Within the core the number densities increase 
from 10'"' to 10^ cm^^. For densities >, 10^ cm~"^, 
however, three-body formation of H2 will become 
the dominant formation mechanism, transform- 
ing all hydrogen into its molecular form (Palla 
etal. 1983). Our chemical reaction network does 
not include this reaction and the solution cannot 
be correct at r <, 0.1 pc. The most interesting ef- 
fect of the three-body reaction is that it will in- 
crease the cooling rate by a factor ^ 10^, lead- 
ing to a further dramatic density enhancement 
within the core. This will decrease the dynamical 
timescales to 4; 100 years, effectively decoupling 
the evolution of the fragment from the evolution 
of its host primordial molecular cloud. Therefore, 
it is a firm conclusion that only the gas within 
these cores can participate in population III star 
formation. 

Omukai & Nishi (1998) have simulated the evo- 
lution of a collapsing, spherically symmetric pri- 
mordial cloud to stellar density including all rele- 
vant physical processes. Coincidentally, their ini- 
tial conditions are very close to our final state. 
Based on their results, we can say that if the cloud 
does not break up, a massive star will be formed. 
Adding a small amount of angular momentum to 
the core does not change this conclusion (Bate 
1998). A third possibility is that the cloud breaks 
up into low mass stars via thermal instability 
in the quasi-hydrostatic phase. Silk (1983) has 
argued that, due to the enhanced cooling from 
the 3-body produced II2 , fragmentation of this 
core might continue until individual fragments are 
opacity limited (i.e. they become opaque to their 
cooling radiation). Exploring which of these sce- 
narios is correct will have to await yet higher res- 
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Fig. 2. — Spherically averaged mass weighted profiles around the baryon density peak shortly before a well 
defined fragment forms (z=19.1). Panel a) shows the baryonic number density, enclosed gas mass in solar 
mass, and the local Bonnor-Ebert mass (w 27 MqTI^^ /^/n). Panel b) plots the molecular hydrogen number 
fraction fn^^ and the free electron number fraction x. The H2 cooling time tn^, the time it takes a sound 
wave to travel to the center, tcross, and the free-fall time i// = [37r/(32Gp)]^/^ are given in panel c). Panel 
d) gives the temperature in Kelvin as a function of radius. The bottom panel gives the local sound speed, 
Cg (solid line with circles), the rms radial velocities of the dark matter (dashed line) and the gas (dashed 
line with asterisks) as well as the rms gas velocity (solid line with square symbols). The vertical dotted line 
indicates the radius at which the gas has reached its minimum temperature allowed by H2 cooling (~ 5pc). 
The virial radius of the 5.6 x 10*' Mq halo is 106 pc. The cell size on the finest grid corresponds to 0.024 pc. 
Note that the simulation box size corresponds to 6.4 proper kpc. 
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olution simulations including the effects of radia- 
tive transfer. It will also be interesting to exam- 
ine the possible effects of molecular HD which, al- 
though much less abundant, is a much more effi- 
cient coolant at low temperatures. 

How many cores are formed in our halo? Be- 
cause our timestep contracts rapidly once the first 
core forms, we are not yet able to answer this ques- 
tion definitively. An earlier, less well resolved sim- 
ulation yielded 5 — 6 cores hy z = 16.5, suggesting 
that multiple cores do form. We speculate that 
the total number of cores to eventually form will 
be proportional to the total amount of cooled gas. 
However, the first star in a given halo will most 
likely always be formed close to its center where 
the dynamical timescale is shortest. The cooling 
timescale at rr^^^ in Figure 2. of ^10^ years 
should roughly correspond to the typical forma- 
tion time of fragments. During this time the prod- 
uct of the first collapsed fragment might already 
be an important source of feedback. Hence, even 
for the question of the efficiency of fragmentation 
it seems that feedback physics have to be included. 

Let us assume the first ~ 200 Mq cores frag- 
ment to form stars with 100% efficiency. If the 
ratio of produced UV photons per solar mass is 
the same as in present day star clusters than 
about 6 X 10^^ UV photons would be liberated 
during the average life time of massive star (~ 
5 X 10^ years). This is about hundred times more 
than the ^ 4 x 10^^ hydrogen atoms within the 
virial radius. However, the average recombination 
time {nkrec)~^ ~ 5 x 10^ years within the virial 
radius is a factor 100 less than the average life- 
time of a massive star. Hence, very small or zero 
UV escape fractions for these objects are plau- 
sible. However, the first supernovac and winds 
from massive stars will substantially change the 
subsequent hydrodynamic and chemical evolution 
as well as the star formation history of these ob- 
jects. A more detailed understanding of the role 
of such local feedback will have to await yet more 
detailed simulations that include the poorly un- 
derstood physics of stellar feedback mechanisms. 

Since the collapsing core evolves on much faster 
timescales than the rest of the halo it seems plau- 
sible that the first star (or star cluster) will have 
a mass of less or the order of the core mass. It 
seems also quite clear that the radiative feedback 
from this star (these stars) will eventually halt fur- 



ther accretion. As a consequence this might sug- 
gest that the formation of very massive objects or 
supermassive black holes is unlikely. This later 
speculation will be tested by yet higher resolution 
simulations we are currently working on. 

Recently Bromm, Coppi, and Larson (1999) 
have studied the fragmentation of the first objects 
in the universe. The results of their simulations 
using a smooth particle hydrodynamics technique 
with isolated boundary conditions disagree with 
the results presented here. Their objects collapse 
to a disk which then fragments quickly to form 
many fragments throughout the rotationally sup- 
ported disk. The efficiently fragmenting disk in 
those simulations originates from the assumed ide- 
alized intial conditions. These authors simulated 
top-hat spheres that intially rotate as solid bod- 
ies on which smaller density fluctuations were im- 
posed. Naturally they find a disk. Also it is clear 
that for a top-hat if the disk breaks up it will do 
so everywhere almost simultaneously. Our results 
with realistic initial conditions do not lead to a 
disk and form the first fragment close to the cen- 
ter of the halo. 

5. Conclusions 

We have reported first results from an ongoing 
project that studies the physics of fragmentation 
and primordial star formation in a cosmological 
context. The results clearly illustrate the advan- 
tages and power of structured adaptive mesh re- 
finement cosmological hydrodynamic methods to 
cover a wide range of mass, length and timescales. 
All findings of AANZ are confirmed in this study. 
Among other things, these are that 1) a signifi- 
cant number fraction of hydrogen molecules is only 
formed in virialized halos at the intersection of fil- 
aments, and 2) only a few percent of the halo gas 
has cooled to T ^ Tyir- 

The improvement of a factor ~ 1000 in res- 
olution over AANZ has given new insights into 
the details of the fragmentation process and con- 
straints on the possible nature of the first struc- 
tures: 1) Only < 1% of the baryons within a 
virialized object can participate in population HI 
star formation. 2) The formation of super massive 
black holes or very massive objects in small halos 
seem very unlikely. 3) Fragmentation via Bonnor- 
Ebert instability yields a ~ 200 M© core within 
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one virializcd object. 4) If the gas were able to 
fragment further through 3-body H2 association 
and/or opacity hmited fragmentation, only a small 
fraction of all baryons in the universe will be con- 
verted into small mass objects. 5) The escape frac- 
tion of UV photons above the Lyman limit should 
initially be small due to the high column densities 
of HI {Nhi ~ 10^^ cm~^) of the parent primordial 
molecular cloud. 6) The first star in the universe 
is most likely born close to the center of its parent 
halo of > 10^ Mq. 
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z=19.1 

3 320pc B-*^ 32pc 




Fig. 3. — Gas density and temperature in the first cosmological objects expected to form in hierarchical 
structure formation scenarios. The upper panels show the log of the baryonic overdensity in a slice through 
the point of highest gas density at a scale of 320 pc (left) and 32 pc (right). The lower panels give the 
corresponding plots of the log of the gas temperature. Additionally the velocity field is also visualized in the 
upper left panel. Note that the computational volume simulated is 20"^ times larger than the left panels. 
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